{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# VARMAX models\n",
    "\n",
    "This is a brief introduction notebook to VARMAX models in statsmodels. The VARMAX model is generically specified as:\n",
    "$$\n",
    "y_t = \\nu + A_1 y_{t-1} + \\dots + A_p y_{t-p} + B x_t + \\epsilon_t +\n",
    "M_1 \\epsilon_{t-1} + \\dots M_q \\epsilon_{t-q}\n",
    "$$\n",
    "\n",
    "where $y_t$ is a $\\text{k_endog} \\times 1$ vector."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import statsmodels.api as sm\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "dta = sm.datasets.webuse('lutkepohl2', 'https://www.stata-press.com/data/r12/')\n",
    "dta.index = dta.qtr\n",
    "endog = dta.loc['1960-04-01':'1978-10-01', ['dln_inv', 'dln_inc', 'dln_consump']]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Model specification\n",
    "\n",
    "The `VARMAX` class in statsmodels allows estimation of VAR, VMA, and VARMA models (through the `order` argument), optionally with a constant term (via the `trend` argument). Exogenous regressors may also be included (as usual in statsmodels, by the `exog` argument), and in this way a time trend may be added. Finally, the class allows measurement error (via the `measurement_error` argument) and allows specifying either a diagonal or unstructured innovation covariance matrix (via the `error_cov_type` argument)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 1: VAR\n",
    "\n",
    "Below is a simple VARX(2) model in two endogenous variables and an exogenous series, but no constant term. Notice that we needed to allow for more iterations than the default (which is `maxiter=50`) in order for the likelihood estimation to converge. This is not unusual in VAR models which have to estimate a large number of parameters, often on a relatively small number of time series: this model, for example, estimates 27 parameters off of 75 observations of 3 variables."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\ProgramData\\Anaconda3\\lib\\site-packages\\statsmodels\\tsa\\base\\tsa_model.py:165: ValueWarning: No frequency information was provided, so inferred frequency QS-OCT will be used.\n",
      "  % freq, ValueWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                             Statespace Model Results                             \n",
      "==================================================================================\n",
      "Dep. Variable:     ['dln_inv', 'dln_inc']   No. Observations:                   75\n",
      "Model:                            VARX(2)   Log Likelihood                 361.038\n",
      "Date:                    Tue, 24 Dec 2019   AIC                           -696.075\n",
      "Time:                            15:05:43   BIC                           -665.948\n",
      "Sample:                        04-01-1960   HQIC                          -684.046\n",
      "                             - 10-01-1978                                         \n",
      "Covariance Type:                      opg                                         \n",
      "===================================================================================\n",
      "Ljung-Box (Q):                61.32, 39.29   Jarque-Bera (JB):          11.38, 2.39\n",
      "Prob(Q):                        0.02, 0.50   Prob(JB):                   0.00, 0.30\n",
      "Heteroskedasticity (H):         0.45, 0.40   Skew:                      0.16, -0.38\n",
      "Prob(H) (two-sided):            0.05, 0.03   Kurtosis:                   4.88, 3.44\n",
      "                            Results for equation dln_inv                            \n",
      "====================================================================================\n",
      "                       coef    std err          z      P>|z|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------------\n",
      "L1.dln_inv          -0.2396      0.093     -2.575      0.010      -0.422      -0.057\n",
      "L1.dln_inc           0.2937      0.449      0.654      0.513      -0.586       1.174\n",
      "L2.dln_inv          -0.1650      0.155     -1.063      0.288      -0.469       0.139\n",
      "L2.dln_inc           0.0834      0.422      0.198      0.843      -0.744       0.911\n",
      "beta.dln_consump     0.9441      0.640      1.476      0.140      -0.310       2.198\n",
      "                            Results for equation dln_inc                            \n",
      "====================================================================================\n",
      "                       coef    std err          z      P>|z|      [0.025      0.975]\n",
      "------------------------------------------------------------------------------------\n",
      "L1.dln_inv           0.0634      0.036      1.775      0.076      -0.007       0.134\n",
      "L1.dln_inc           0.0835      0.107      0.779      0.436      -0.126       0.293\n",
      "L2.dln_inv           0.0105      0.033      0.318      0.750      -0.054       0.075\n",
      "L2.dln_inc           0.0365      0.134      0.272      0.786      -0.227       0.300\n",
      "beta.dln_consump     0.7688      0.112      6.848      0.000       0.549       0.989\n",
      "                                  Error covariance matrix                                   \n",
      "============================================================================================\n",
      "                               coef    std err          z      P>|z|      [0.025      0.975]\n",
      "--------------------------------------------------------------------------------------------\n",
      "sqrt.var.dln_inv             0.0434      0.004     12.284      0.000       0.036       0.050\n",
      "sqrt.cov.dln_inv.dln_inc  5.848e-05      0.002      0.029      0.977      -0.004       0.004\n",
      "sqrt.var.dln_inc             0.0109      0.001     11.226      0.000       0.009       0.013\n",
      "============================================================================================\n",
      "\n",
      "Warnings:\n",
      "[1] Covariance matrix calculated using the outer product of gradients (complex-step).\n"
     ]
    }
   ],
   "source": [
    "exog = endog['dln_consump']\n",
    "mod = sm.tsa.VARMAX(endog[['dln_inv', 'dln_inc']], order=(2,0), trend='n', exog=exog)\n",
    "res = mod.fit(maxiter=1000, disp=False)\n",
    "print(res.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From the estimated VAR model, we can plot the impulse response functions of the endogenous variables."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAwcAAADcCAYAAAAlSd0+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xl8XPV57/HvMyPJsuVdli28yJbxxhY7RDYYy7EIECCXbG4MKQTC4rhNmqS0N83S0tbNJdDb5ibktklvXZsEAklKQgLZICQhXgQ2wTaQsHnBC8jYWJaN5d2S5rl/nCN5JI2k0Wg5I+nzfqGXzpz5nd95ZjzGv2d+v3Mec3cBAAAAQCzqAAAAAABkB5IDAAAAAJJIDgAAAACESA4AAAAASCI5AAAAABAiOQAAAAAgieQAADBAmdl1ZrYk6jgyZWZXmtknoo4DQP9i1DkAAAxEZjZD0hpJc929Kup4OsvMxkjaKOkj7r4x6ngA9A/MHADok8xsuZkdNbP9ZrbHzP5n1DENJGb2HTO7OVv6yYS7b5X0Z5KWp3rezCrMbHWLfbvMbEpXz21mz5rZpK704e4HJH1E0pe7Gg8ANCI5ANCX/bu7j5W0QNIXzGx21AH1JDObEtVAOpuY2Rwz+1B39OXuP3X3pd3RVyfPO9fd3+iGfja6+/u6IyYAkEgOAPQD7r5L0gZJMyMOpadNkXRzxDFkgzmSuiU5AAA0R3IAoM8zsxJJZZK2hI9vMbMdZrY3+YJNM/snM3sz3P+pcN9yM/uZmb1iZq+b2fuT2v9NuG+LmV2d1P7rZva4mdWY2TfC/TEzWxH2/bqZfTipn78N9+1u7L+99m28xnWSfizpEjPbZ2b/ldTP/wmXVr1gZnM76GeImf0kPO+LZvbODtoPDtvvM7OtZrYg6emJZva0mR00s88l9f+dsP9KM5sW7h9jZo+Gy8CeMrPSFueJm9nPzez2DuLZJekbkq4LY/qH9s7bWWY2LVzys0vS9Wm0/46Z3dHyfUjjuGbLk9rqx8x+aGbXhNuDwtc3KJPXBgDpIDkA0Jd92sz2S9om6V/d/QUzO0/SXytIFuZIWm5m48xstKQvSTpP0vmSLkvq5zxJl0h6v6R7w4Hm5ZI+Lmm2pA9L+raZjQvb3yrp7yW9Q9Kfh33PkXSNpMmSrpR0hSSFScVlks6R9F5JK8wst632bXH3hZIWS3ra3YvdvTHpuTXs62xJfyXphx0MHq+UtFfSeEl3hu9Je66WNCFsv0zS5UnP/ZmkG8M+/z7c97eS4uExKyV9L9z/DUkvSBon6eeSvtLiPN+U9LK739NeMO4+RdJfSvrv8H1oXG/f1nk76/+Gx04NY01HqvchE6n6+ZGCPwNJWiTpSXc/1YVzAEC7SA4A9GX/rmCpzVFJvwj3vUfBwO5lBYPRoQqWGx1WMLPwdQUD8Y8n9fOIux9y9xckvSVphoIB2QPh/pclPSNpYdj+Z+7+rLvvkbRP0nBJr0lKSPrX8PjGb8AvlzQ3fH6NpAIFA+222nfW1ZL+y91PuvuT4eu8oK3G7v4TSY9K+qqkf5A0toP+X5A0SdJdkvIk/a+k5+5z99cU3DFneFI833T3hLt/R9IMC+6qc1W43939bndP/lb+Mwq+pV+fzgtuQ1vn7axLFPy5JyR9O81jUr0PmUjVz88lVYTbVytIFgCgx5AcAOjT3P24pHslfSrcZZLuD79VLpY0UdIGd29QMEj/kYJvYJ8zs7ykYxrFFAzaJSn5Xs+e9Pi1Fvvl7oclnStpnYKB7q+T+v5KUjwlkva00z4TbcXZipn9naTPK0hU/rbDjoPB6jskvSjpc2o+YH4tbNPyfO0+NrORZnZj0q4aBYnXnWYW7yim9sLt4HE6LOm4RHsNk7T1PnRWq37c/ZikrWY2U8F79HgXzwEA7SI5ANAf/LukG82sQNKTkq42s2IzG6bgm+9zLbin/W/Dny9IKpZUGB7/ITMbFa6/L1SwTOkxSTeEA9lZki6SVBm2bzUINLPLFCQpjyhYqjPPzEzSbyRda2bDzaxxxmBkO+3bc0DBOv94GG88jPO2cD36IkkjFQzk2zJf0g/DuNq9ziF8XbdIukPBUpt/lnRx0tOpBsOPSfpkeC3EjZK2unuNpF9J+mTY5k8V3IKz0ffCWZvNkm7rKCYF78PkML7G2YG2zttZv5f0p2YWk3RTmsd0V8Ggtvr5kYLkd7u7n+imcwFASiQHAPo8d98taa2kG9z9RQVLX9YrWFr0TXd/Pryn/TpJOyVtVXAb1L1hF88pGPg/KulWdz/h7r+R9F1Jf5D003D/W+2EsUbSEUl7wvN8PlxC80sFCcCLkp6S9Jnw/vQp23fwOl9UMKjfE/aXpyDB+KOkHQrWyy/pYE36txQsJ3pZ0tuSpnbwbf2PJJUqWD51r4JEoT13Kfj2fY+kP9eZi3pvl/QuM9sr6aOSPpvi2OWS7giTvPb8StJhM3tLwfvR3nk767MKruN4Q1K2rO3/maSlkh6OOhAA/R8VkgEMaGa2XJLcfXm0kQAAED1mDgAA6EZm1ni72ZY//y/q2ACgI8wcAAAAAJDEzAEAAACAEMkBAAAAAElSTtQBtGfMmDE+ZcqUqMMAAAAA+qxNmzYdcPeidNpmdXIwZcoUbdy4MeowAAAAgD7LzHan25ZlRQAAAAAkkRwAAAAACJEcAAAAAJCU5dccAAAAAC3V1dWpqqpKJ0+ejDqUrJKfn6+JEycqNzc34z5IDtpx//pdOnSsTn95+fSoQwEAAECoqqpKw4YN05QpU2RmUYeTFdxdNTU1qqqqUmlpacb9sKyoHS/uOaz/WLNdh46djjoUAAAAhE6ePKnCwkISgyRmpsLCwi7PppActOO28qk6WZfQ937/etShAAAAIAmJQWvd8Z6QHLRjZvEwLZw+Rt95epdO1TdEHQ4AAADQo0gOOrB04VRVHzmln7+wN+pQAAAAgB7FBckdePf0MZo+dqhWVe7U4gsnMIUFAACQRf7pZy/p5Tdru7XPc8cP1z++/7y02y9fvlwVFRVNP6tXr87ovLfffrvuueeejI7tLswcdMDMtHRhqV7eW6v1O2qiDgcAAAD9VNSJgcTMQVo+OGeC/uXxLVq1bqcuOXtM1OEAAAAg1Jlv+LvToUOHtGTJEjU0NMjdVVFR0arN8uXLVVdXp3Xr1qm2tlaPP/64iouL2+wzedYh1bGrVq3Seeedpw996EO6++67NW3aNC1ZsqRbX1dGMwdmtsrM1pvZHZ1tY2bjzOy5TM4blfzcuG6cP1m/fXW/Xqs+GnU4AAAAiNiKFSt0zTXX6He/+127Rce2b9+utWvXavHixXryySc7dY6Wxy5ZskSPPfaYJGnt2rV63/ve16XXkEqnkwMzWywp7u7zJU01s1YVwjpo81VJgzMNOCofu3iy8nJiurdyZ9ShAAAAIGI7d+7U7NmzJUllZWVttrvpppskSSUlJTp9unO1s1oeO2PGDFVVVam2tlYjR45UQUFBhtG3LZOZgwpJD4XbT0gqT7eNmb1H0jFJ+9rq3MyWmdlGM9tYXV2dQXg9Y8zQQfrwnAl6eHOVDlIUDQAAYEArKSnRSy+9JEl6/vnn22zXlQF8qmPnzZune+65Rx/4wAcy7rc9mSQHBZL2hNsHJY1Lp42Z5Un6e0lfbK9zd1/h7mXuXlZUVJRBeD3ntoWlQVG0Z3ZHHQoAAAAitGzZMj388MOqqKhQbW333i2pPUuWLNE999yja665pkf6N3fv3AFm35D0fXffEC4fmuXud3XURsHFz6+4+w/NbLW7V3R0rrKyMt+4cWOn4utpN937e72yt1aVX7hUg3LiUYcDAAAw4Lzyyis655xzog4jK6V6b8xsk7u3vfYpSSYzB5t0ZinRbEm70mxzuaS/MLPVkuaY2coMzh25peWlFEUDAABARpLrIVRUVOiDH/xg1CE1k8mtTB+RtM7Mxku6WtJHzexOd7+jnTYXu/v3Gp8MZw6WdiXwqCycPkYzxg3VSoqiAQAAoJMyLZDWWzo9c+DutQouON4g6VJ3f6FFYpCqzeEWz1dkGG/kzExLy6fqlb21Wv8aRdEAAADQf2RU58DdD7n7Q+7e5l2H0mnTV31gzniNGZqnldzWFAAAAP1IRsnBQJefG9fHLp6sJ1/dr+37KYoGAACA/oHkIENNRdGeYvYAAABgIFu+fHnTtQQVFRUZ93P77bd3T0BdkMkFyVBQFG3xOyfo4U1V+tx7Z2p0QV7UIQEAAAw8j31R2vfH7u2z+ALp6n/u3j7TcM899/T6OVsiOeiCW8tL9YNn39D3ntmtT79netThAAAAoJccOnRIS5YsUUNDg9w95YzB8uXLVVdXp3Xr1qm2tlaPP/64iouL2+yzoqKiaQYi1bEjR47UzTffrKqqKo0cOVIPPfSQhgwZ0q2vi+SgC2aMG6ZFM4p03/rd+sS7p1IUDQAAoLdF8A2/JK1YsULXXHONbr/9dl1xxRVtttu+fbvWrl2rL3/5y3ryySd1/fXXp32OlsceOHBAs2fP1g9+8AN9+9vf1osvvqh58+Z1x8tpwjUHXbR0YVAU7WcURQMAABgwdu7cqdmzZ0uSysraLj580003SZJKSkp0+vTpTp2j5bGvvvpqUzJw8803a+7cuZmE3i6Sgy4qnzZGM8cN08p1O+TuUYcDAACAXlBSUqKXXnpJkvT888+32a6goCDjc7Q8dtasWXr22WclSXfddZdWrlyZcd9tITnoIjPTbeWlenXfET1NUTQAAIABYdmyZXr44YdVUVGh2traXjnnJz7xCW3evFkVFRXavHmzbrzxxm4/h2Xzt91lZWW+cePGqMPo0Mm6BpX/7yd1wYQR+vYt3bvuCwAAAM298sorOuecc6IOIyulem/MbJO7t732KQkXJHeD/Ny4brx4ir7+m63avv+Ipo0dFnVIAAAAyEIt72o0YsQIPfroo9EEkwLJQTf52MUl+ubq7VpVuUt3L74g6nAAAAD6NXeXmUUdRqc13qq0J3THiiCuOegmhUMH6U8unKAfb67SwWOduxIdAAAA6cvPz1dNTQ03g0ni7qqpqVF+fn6X+mHmoBvduqBU3//9G3pww2595jKKogEAAPSEiRMnqqqqStXV1VGHklXy8/M1ceLELvVBctCNpo8bpoqZQVG0ZYsoigYAANATcnNzVVpaGnUY/RLLirrZ0vKpOnD0lH76/JtRhwIAAAB0CslBN1swrVCziodpVeVO1sEBAACgTyE56GZmplvDomhPbacoGgAAAPoOkoMe8ME54zVm6CCtrNwRdSgAAABA2no1OTCz0WZ2hZmN6c3z9rZBOXHdNH+yVm+p1vb9R6IOBwAAAEhLRsmBma0ys/Vmdke6bcxslKSfS5on6XdmVpRRxH3EDReVaFBOTKsqd0UdCgAAAJCWTicHZrZYUtzd50uaamatbujfRpt3SPprd/+KpF9JurBroWe3wqGDtPjCifrx5irVHD0VdTgAAABAhzKZOaiQ9FC4/YSk8nTauPsad99gZu9WMHuwPlXnZrbMzDaa2ca+XtjitvIpOlWf0IPPvB51KAAAAECHMkkOCiTtCbcPShqXbhszM0nXSTokqS5V5+6+wt3L3L2sqKhvrzyaNnaYLp1ZpPvX79LJuoaowwEAAADalUlycFTS4HB7aBt9pGzjgb+Q9AdJH8jg3H3ObeVTdeDoaf30BYqiAQAAILtlkhxs0pmlRLMl7UqnjZl9wcxuCveNlPR2Bufuc5qKoq2jKBoAAACyWybJwSOSbjSzr0m6VtJLZnZnB21+IWlFuG+tpLiCaxH6PTPTbeWl2vLWEVVuPxB1OAAAAECbOp0cuHutgguON0i61N1fcPc7Omhz2N0PufsV7v5ud/+UD6Cv0T/QWBRt3c6oQwEAAADalFGdg3Cg/5C77+tKm4FiUE5cH58/WWu2VmvbWxRFAwAAQHbq1QrJA9kNF0/WoJyY7n2K2QMAAABkJ5KDXjK6IE9/8q6JenjzHoqiAQAAICuRHPSiWxeU6nR9Qg9soCgaAAAAsg/JQS+aNnaoLp1ZpO9uoCgaAAAAsg/JQS9bujAsivY8RdEAAACQXUgOetklZwdF0VZW7qAoGgAAALIKyUEvMzMtXThVW986SlE0AAAAZBWSgwi8f/ZZKhpGUTQAAABkF5KDCCQXRdtKUTQAAABkCZKDiFx/0WTl58Z0byWzBwAAAMgOJAcRGV2Qpz+5cKJ+/NweHaAoGgAAALIAyUGEbi1vLIq2O+pQAAAAAJKDKJ1dNFTvmTVW312/m6JoAAAAiBzJQcSWlpeq5thpPfr8nqhDAQAAwABHchCx+WcX6pyzhmtV5U6KogEAACBSJAcRMzMtLS/V1reOat02iqIBAAAgOiQHWeD9s8dr7LBBWsltTQEAABAhkoMskJcT08cvmaK1W6u1ZR9F0QAAABANkoMscf28EoqiAQAAIFIZJQdmtsrM1pvZHem2MbMRZvaYmT1hZj8xs7xMg+6PRoVF0X7y/B5VH6EoGgAAAHpfp5MDM1ssKe7u8yVNNbPpaba5QdLX3P29kvZJuqprofc/FEUDAABAlDKZOaiQ9FC4/YSk8nTauPu33P3X4b4iSftTdW5my8xso5ltrK6uziC8vuvsoqG6bNZYPbCBomgAAADofZkkBwWSGit2HZQ0rjNtzGy+pFHuviFV5+6+wt3L3L2sqKgog/D6ttsWUhQNAAAA0cgkOTgqaXC4PbSNPlK2MbPRkv5N0q0ZnHdAmD+1UOeeNVwr11EUDQAAAL0rk+Rgk84sJZotaVc6bcILkH8o6UvuzqL6NpiZli4s1bb9R7WWomgAAADoRZkkB49IutHMvibpWkkvmdmdHbT5haTbJF0o6e/MbLWZXdeFuPu1a94RFkVbtyPqUAAAADCAdDo5cPdaBRccb5B0qbu/4O53dNDmsLv/h7uPcveK8Oe/ux5+/9RYFG3dtgMURQMAAECvyajOgbsfcveH3H1fV9qgbTdcFBRFW1XJ7AEAAAB6BxWSs9TIIXn6yLsm6pHn3qQoGgAAAHoFyUEWu3VBqU43JPRdiqIBAACgF5AcZLGpRUN1+Tlj9SBF0QAAANALSA6y3G3lU1Vz7LQeeY6iaAAAAOhZJAdZ7uKpo3Xe+OFaWUlRNAAAAPQskoMs11gUbfv+o1qztTrqcAAAANCPkRz0Af/jgvEaN3yQVlXujDoUAAAA9GMkB31AXk5MN80PiqK9uq826nAAAADQT5Ec9BE3XFSiwblxrVrH7AEAAAB6BslBH9FYFO3R5ymKBgAAgJ5BctCH3LJgiuoSFEUDAABAzyA56EOmFg3VZbPG6QGKogEAAKAHkBz0MUsXlurgsdP6CUXRAAAA0M1IDvqYi0pH6/wJw7WqcqcSCYqiAQAAoPuQHPQxZqbbysOiaNsoigYAAIDuQ3LQBzUVReO2pgAAAOhGJAd9UF5OTB+/ZIoqtx/QK3spigYAAIDuQXLQR10/LyiKdm8lswcAAADoHhknB2a2yszWm9kdnWljZuPMbF2m50Vg5JA8LSkLiqLtP3Iy6nAAAADQD2SUHJjZYklxd58vaaqZTU+njZmNknSfpIKuBI3ALQtKVZdI6IH1FEUDAABA12U6c1Ah6aFw+wlJ5Wm2aZB0naQ2F8qb2TIz22hmG6uruRtPe0rHFOjyc8bpuxRFAwAAQDfINDkokNRYheugpHHptHH3Wnc/3F7H7r7C3cvcvayoqCjD8AaOpeWlOnS8Tj/eTFE0AAAAdE2mycFRSYPD7aFt9JNOG3TRvKaiaDsoigYAAIAuyXTAvklnlhLNlrQrwzboIjPT0vKpeq36mNZsZRkWAAAAMpdpcvCIpBvN7GuSrpX0kpnd2UGbX2QeJtrzvgvOUvHwfK2s3BF1KAAAAOjDMkoO3L1WwQXHGyRd6u4vuPsdHbQ5nPRcRYbxIoXGomhPba+hKBoAAAAylvF1AO5+yN0fcvd9XWmD7tFYFG0VRdEAAACQIS4S7idGDMnVtWUT9ejze7S/lqJoAAAA6DySg37klgWlqk+4vruBomgAAADoPJKDfmTKmAJdcc44PbBht06cpigaAAAAOofkoJ+5rbEo2nNVUYcCAACAPobkoJ+ZVzpaF0wYoVWVOymKBgAAgE4hOehnzExLF5ZqR/Uxrd66P+pwAAAA0IeQHPRDjUXRuK0pAAAAOoPkoB/Kjcd084KgKNrLb1IUDQAAAOkhOein/nRuiYbkURQNAAAA6SM56KeComiT9NMXKIoGAACA9JAc9GO3LJii+oTr/vUURQMAAEDHSA76scmFYVG0ZyiKBgAAgI7lRB0AetbShVP1xMtv6eHNVfrYxZM7PsBdOlYt1e6RDu8Jf4cF1UZMkkZMkIZPkEZMlIYUSmY9+wKQFU7WNejZXQe1eku1KrcdUF0ioYmjhmjiqMGaFP4OfoZozNA8GZ8LAAD6JJKDfm7ulFF6x8QRurdyp66fO0mxU28nDfyrkhKA8HHtm1LD6eadxPOC3y335+SHicIEafjE4PeIiWe2h0+Q8of3zgtFt3J37ao5rjVb9mvN1mqt31Gjk3UJ5eXEdFHpaA0dlKOqQyf0x6q3deh4XbNj83NjmjBycFPyMHHUEE0afeZxYQHJAwAA2YrkoD85Wdtq4G+1e7TKXlNt7W75XW9LDSeaHxPLkYaNDwbzE8qkc8aHA/ykQX/BmGBG4fiBYBahcTahaXuPtHONdGSv5Inm/Q8a3rq/piQiTCBy83vvPUKbjp2q1/rXarRma7XWbK3W6wePS5JKxxToo3NLtGhmkS4uLdTgvHiz446eqteeQydUdei4qpr9PqE/tJE8TGwx25A8AzGa5AEAgMiYu0cdQ5vKysp848aNUYeRHU4fPzMob+ub/9NHmh9jMWlosRLDJ+h3e3N1anCx3ldedmZZ0PAJ0tCxUiye+pyd1VAfJAipkofGWI8faH3ckDFhwjApKYkIYxwxURpaLMXJY7ubu2vLW0e0ZkuQDDy766DqGlxD8uK65OxCLZpRpHfPKNLkwoIunefIyTrtefuEqg42Txyq3g62326RPAzOjbdKHJoSiNFDNGpILskDAACdYGab3L0snbaMuLJB/anWa/xbPj75duvjCsYGA+nCadLUCmn4+OYD/2HFUjxXMUnb17ymux97Vb+YXK7zxo/omdcRz5FGTgp+2lJ3Ili6lPwaD78RbNe8Ju1cK51qUbjNYtKws5IShgnNly6NmBTMbjBg7NDh43Wq3H5Aa7YGy4Xeqj0lSZpVPEy3LijVohlFeteUURqU000Jo6Rh+bmaVZyrWcWpl5gdOVl3JmFoMfuw+fW3dfhE8+RhSF68ReLQPIkgeQAAIHPMHPS0hrrg2/SUA/9wjf+x6tbHDR7dYhA8vvmAePh4KWdQ2mEcPlGn+Xf/VledX6yvXTunG19gDzh5OMX71WImouFU82Pig4L3pGkJ08TWy5jyeygpymKJhOvFNw9rzZZqrd5aredeP6SES8Pzc7RwelHT7EDxiOxd2lV7si5ctnRCbxxsuXTpuGpP1jdrX5AXbzNxmDhqsEaSPAAABhhmDnpLokE6+laKJT7hILb2TenIPkktErBBI84M8se/88za+8bB7PDxUt6Qbg11xOCgKNqDz+zWF66apXHDs3cwqPwRwc+4c1M/7y4drwlmHFIlEbsqw+sfWty+NW9Y0mxDy2VMYRKRO7jnX18PO3D0lNZtq9aaLdVau+2ADh47LTPpHRNG6NOXTtOimUWaPXGkcuJ9407Gw/NzNfysXJ1zVuqZh8Mn6lpc83BCb4Tbv995UEdOpU4eki+STk4eRgwmeQAADFwZzRyY2SpJ50r6hbvfmW6bdI5LFvnMwbEa6e3dKdb4vxnsO7JXSjQfeCi3oPkAtOXAf8QEadCwSF7O7ppjqvjqan2q4mz9zZWzIomh1zTUh4lbVevErTGRSDVjM6Sw+dKslklE/gjJ4sF1Gk2/Y5EuaapvSOi5N95uunbgj3sOS5IKC/L07hnB7MDC6WNUODT9mab+5PCJuiBxOHhMew4e1d5Dwc++g0f11uFjOnnqlHKUUNwalKMGDc8zjR+eq7OG5eisYbkqHpajcQW5Gjs0rrEFORqS4zJvCP7uu4d/9pbGb4W/Y504JkUfFsvw2PAz2qXzt9OHlPT3oOXjVPs6+bgX/o41JFx1DQmdbkiorj6huoakxw0J1dW7Tjc06HR9sL/x53SDq64+IQ+jtvBtiYUxm9mZ/TLFmt5KC/edaROLBW3UYn9jf00fo6RjY7EzbVr2GUtqq+QYYs37SD421iLexmOTz5N8rJodY41/ck3t1Ozcjdsp2jXFkfQ7OfbwuDPbZx4DaF+PzhyY2WJJcXefb2b3mtl0d9/WURtJF3R0XNb58VLptSfPPI4POjPgn1Leeo3/iAlS/she+UcsE5MLC/Tec8fpwWde16cvnd7qrjP9SjwnHNhPkHRR6jZ1J6Ujb4bXPbRIIg7tlnY9JZ06nOYJrUXCEA//lW+ZSLS1P5aiXVv7YzpRL+0/Wq99R+q07+hpnao3TbKYPjN0sIpLCzR+VIEKhw2WxeLSgbhU05XzpdjvHgyOE/XBDFrTdlv7Ovu4O/pokLxBIxL1GpGo13kt76QV/rEp1STakfAHWcvDwWTL36mf92bPu6xpPtcluVvTtiTFZBokKU/N96vF47bOHbT0Nh+3bt2ybXr9dKVt63N2V7+9x9VqXr57Zec/5W2/5i6+Gak+x92hp/rta7a8936dv+CaqMNISybLiiokPRRuPyGpXFLLQX6qNu9M4ziZ2TJJyySppKQkg/C60YLbpblL+1XRr6ULp+pXL3WiKFp/lpsvjZ4a/LTl1JGkWaMq6dTRYLlSOPBUItHicUNwO9dmj5P2p3yuvT5OS4kGJRINOn7ytI6dPKUTp06rvr5eMSU0IeaakWvKHywNirtinpDebpAOtdFvz/5T2rFYbnD73FhOkGQ0baf5OGdQO8/nBIlXp/psv41bXMfrTfuP1Qc/Rxu0L0zK9tbW6Xh9g5SQGhINSrgkTyiRcHn4O+EJubs8ESQmFg5LY+E/l2ZnhqmNz1mKbbXYH0sayxqFAAAL80lEQVR6vuVzLbeVYn+s1XmluLniJsVMrX7HzJVjJktqk2ONx0kJD15rQ0JKJBJKuKshkVBDwtWQUJBMKnlQ2fhYzfZbW/ut+f7G43NiUtxMcZPiMVNOLPgGOycWfDMet6THYZvGtsmPY2EfsabnlXJ/zJKP86bHsZgFfyYWk8sV/pf0CvzMn2Kz57xZyuDhS/NWyUz4u8WxjU80S3aS/oo3ntP9TC9NnwpP6jfFeYLfyfE3HuctXpuS4m59/qbtFqsUPGmjWbsWO9Jt16y9J7VpGUeK4Jr/mXXQrun1tHo66bWn+v9s85mR1M+22JdiZ8p2qbvs+PikWZieitmaZlcxafzZUYeQtkySgwJJe8Ltg5IuTLNNOsfJ3VdIWiEFy4oyiK/7TF0U6el7QtnkUZrdWBRtXoliMf7StmvQMGnsrOCnl+2uORbUHNhSradfq9GJugblxWOaVzpai2YUadHMIpWOHdq5aXX3thOXzuy3TAbhfeMah2Sm4H9cpeFPVyQSrgb3cMAcbCcSrvqENz1X3+DhwDr4XR+2TSSk+qZBt1r10dDYT3L/Lc/R3nNNsSkc1Ev17jrV6vnUx8ZjptycmPLiMeXGTbnxmHLjMQ3KiTVt5+ZY+Hzjjymv6ZiYcnPCfUlt8nIs/J20Lzw2HjOWlABAD8gkOTgqqfGqzaGSUv2Ln6pNOsehh5mZbls4VZ/9/nP63Zb9uuyccVGHhNDx0/XasKOm6dqBXTVBEbLJhUN0bdnEoAjZ1EINyevCfQTMwpoR3IugtwXfLpty+/FqPgBA35fJCGGTgiVBGyTNlrQlzTZVaRyHXnD1+cUaPyJfK9ftJDmIkLtr2/6jTcnA73ce1OmGhAbnxjX/7ELdEtYdmDKma0XIAAAA0pVJcvCIpHVmNl7S1ZI+amZ3uvsd7bS5WMGCtpb7EIHceEw3L5iiu375ql7cc1jnTxh49/+PyuETdXp6+4FgudDWau09fFKSNGPcUH38kslaNGOsyqaMUj5fLwMAgAhkeivTUZKukLTW3fel2yad45JFfivTfuzwiTpdcvdvdeV5xfradVleFK0PSyRcL++t1Zqt1Vq9Zb82v/62GhKuYYNyVD59TFMRsvEj+359BQAAkJ16vAiaux/SmTsPpd0mnePQO0YMztWSskl6YMNuff6qWVldIbevqTl6SpXbD4RFyKp14OhpSdL5E4brk4vO1qKZRZozaaRy+0gRMgAAMHBwVeIAduuCUt23fpfuX79Ln7+qnxdF60H1DQm9UHWmCNkf9hyWuzRqSG5SEbIiFQ0bmEXIAABA30FyMICVFA7RlecW63u/f12ffs+0rt0FZ4DZd/ik1obXDazbVq3ak/WKmfTOklH6q8tnaNGMIp0/YYTi3CoWAAD0IYwGB7ilC0v1+Ev79PDmPbpxoBdFa8ep+gZt2nVIa7YFdQde3ReUzx03fJCuOr9Yi2aMVfm0MRoxJDfiSAEAADJHcjDAvWvyKM2eNFL3Vu7UDRRFa9J4m9G1W6tVuf2ANuyo0cm6hHLjprlTRutLV8/SoplFmjluGIWYAABAv0FyMMCZmZaWl+oz339OT766X5efO3DrHjReSLxu2wGt21att2pPSZKmFhXoo3NLtHD6GF00tVBDB/HXBgAA9E+McqCrzy/WhJGDtbJyx4BKDhqXCq0Nk4GX3qyVJI0ckqsF08bo3dPHqHx6kSZwm1EAADBAkBxAOfGYbr5kir7yy1f6dVG0xqVCjTMDjUuFcmKmCyeP0ufeO0MLp3MhMQAAGLhIDiBJum7eJN3zm61aVblTX+9HRdFYKgQAAJA+RkSQJA3Pz9W1cyfpu+t36wt9uCgaS4UAAAAyR3KAJrdcUqr7nu5bRdHaWyr0LpYKAQAAdArJAZqUFA7RlecV68FnsrsoGkuFAAAAegajJzSzdGGpHntxnx7eVKUb50+JOhxJLBUCAADoLSQHaObCklGaM2mkVlXu1A0XTY6kKFrLpULP7DioE3UNLBUCAADoYSQHaMbMtHRhqT79vef021f364peqnvQ3lKh6+ZOYqkQAABAL2CkhVauOi8sirZuR48lB8lLhSq3V+vFPSwVAgAAiBrJAVpJLor2x6rDumBi14uidbRU6G+unKnyaWNYKgQAABAhkgOkdKYo2g7d89F3ZtQHS4UAAAD6FkZlSGl4fq6um1ui+9fv0hevPietomgsFQIAAOjbOp0cmNkqSedK+oW739mZdmY2TtKP3H1hhvGiF92yYIq+8/RO3bd+l76Qoiiau2v7/qNNtxhNtVRo4fQxOm88S4UAAAD6gk4lB2a2WFLc3eeb2b1mNt3dt6XTTtIBSfdJKuiWyNHjJo0eoqvOL9aDG3br05dOU8GgnDaXCp3NUiEAAIA+r7MjuApJD4XbT0gql9QqOWij3cOSrpP0aHsnMLNlkpZJUklJSSfDQ3e7rXyqfvnHffrM95/T/iMnWSoEAADQj7WbHJjZf0qambRrkaRV4fZBSRe2cWiBpD3J7dy9Nuyz3YDcfYWkFZJUVlbm7TZGj3vX5FGaVzpaa7dWs1QIAACgn2s3OXD3P0t+bGbfkNT4FfFQSbE2Dj2aZjv0AfffOk8Jdw3JY6kQAABAf9bZQfsmBUuEJGm2pF1dbIc+ID83TmIAAAAwAHR2xPeIpHVmNl7S1ZIuNrNzJV3v7ne0165bogUAAADQYzo1cxBeN1AhaYOkS939sLu/3CIxSNku6bmKLsYMAAAAoAd0eq2Iux/SmTsRdbkdAAAAgOzAhcIAAAAAJJEcAAAAAAiZe/aWEjCzakm7Iw5jjILqzoDE5wHN8XlAMj4PSMbnAcmi/jxMdveidBpmdXKQDcxso7uXRR0HsgOfByTj84BkfB6QjM8DkvWlzwPLigAAAABIIjkAAAAAECI56NiKqANAVuHzgGR8HpCMzwOS8XlAsj7zeeCaAwAAAACSmDkAAAAAECI5AAAAACCJ5KBdZrbKzNab2R1Rx4JomdkIM3vMzJ4ws5+YWV7UMSFaZjbOzJ6LOg5kDzP7lpm9P+o4EB0zG2VmvzSzjWb2n1HHg+iE/0asC7dzzexnZvaUmd0adWwdITlog5ktlhR39/mSpprZ9KhjQqRukPQ1d3+vpH2Sroo4HkTvq5IGRx0EsoOZLZRU7O4/izoWROpGSQ+G97MfZmZ94r726F5mNkrSfZIKwl2fkbTJ3RdI+oiZDYssuDSQHLStQtJD4fYTksqjCwVRc/dvufuvw4dFkvZHGQ+iZWbvkXRMQaKIAc7MciX9l6RdZvbBqONBpGoknW9mIyVNkvRGxPEgGg2SrpNUGz6u0Jkx5VpJWZ00khy0rUDSnnD7oKRxEcaCLGFm8yWNcvcNUceCaIRLyv5e0hejjgVZ4yZJL0v6F0nzzOwzEceD6FRKmizps5JeUTB+wADj7rXufjhpV58aU5IctO2oziwZGCreqwHPzEZL+jdJWb9eED3qi5K+5e5vRx0IssY7Ja1w932SHpB0acTxIDr/KOnP3f3Lkl6VdEvE8SA79KkxZVYHF7FNOrOUaLakXdGFgqiF3xb/UNKX3H131PEgUpdL+gszWy1pjpmtjDgeRG+7pKnhdpkk/h8xcI2SdIGZxSVdJIliUpD62JiSImhtMLPhktZJ+q2kqyVd3GKKCAOImX1S0l2SXgh3/Ye7/3eEISELmNlqd6+IOg5EK7y48F4FSwVyJX3E3fe0fxT6IzObJ+nbCpYWrZf0YXc/Gm1UiErjvxFmNlnSLyX9RtIlCsaUDdFG1zaSg3aEV5tfIWltOF0MAAAAdIqZjVcwe/CrbP+ymeQAAAAAgCSuOQAAAAAQIjkAAAAAIInkAAAAAECI5AAAkDEzm2Nmc6KOAwDQPUgOAABdMSf8AQD0A9ytCACQETO7W9KHw4d73P2yKOMBAHQdyQEAIGNmdrMkuft3oo0EANAdWFYEAAAAQBLJAQCga05IGiJJZmYRxwIA6CKSAwBAV/xa0mIze0rSwqiDAQB0DdccAAAAAJDEzAEAAACAEMkBAAAAAEkkBwAAAABCJAcAAAAAJJEcAAAAAAiRHAAAAACQJP1/yVGOhOjWUt8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 936x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "ax = res.impulse_responses(10, orthogonalized=True).plot(figsize=(13,3))\n",
    "ax.set(xlabel='t', title='Responses to a shock to `dln_inv`');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Example 2: VMA\n",
    "\n",
    "A vector moving average model can also be formulated. Below we show a VMA(2) on the same data, but where the innovations to the process are uncorrelated. In this example we leave out the exogenous regressor but now include the constant term."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\ProgramData\\Anaconda3\\lib\\site-packages\\statsmodels\\tsa\\base\\tsa_model.py:165: ValueWarning: No frequency information was provided, so inferred frequency QS-OCT will be used.\n",
      "  % freq, ValueWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                             Statespace Model Results                             \n",
      "==================================================================================\n",
      "Dep. Variable:     ['dln_inv', 'dln_inc']   No. Observations:                   75\n",
      "Model:                             VMA(2)   Log Likelihood                 353.887\n",
      "                              + intercept   AIC                           -683.775\n",
      "Date:                    Tue, 24 Dec 2019   BIC                           -655.965\n",
      "Time:                            15:05:52   HQIC                          -672.671\n",
      "Sample:                        04-01-1960                                         \n",
      "                             - 10-01-1978                                         \n",
      "Covariance Type:                      opg                                         \n",
      "===================================================================================\n",
      "Ljung-Box (Q):                68.79, 39.25   Jarque-Bera (JB):         12.51, 13.94\n",
      "Prob(Q):                        0.00, 0.50   Prob(JB):                   0.00, 0.00\n",
      "Heteroskedasticity (H):         0.44, 0.81   Skew:                      0.05, -0.49\n",
      "Prob(H) (two-sided):            0.04, 0.59   Kurtosis:                   5.00, 4.87\n",
      "                           Results for equation dln_inv                          \n",
      "=================================================================================\n",
      "                    coef    std err          z      P>|z|      [0.025      0.975]\n",
      "---------------------------------------------------------------------------------\n",
      "intercept         0.0182      0.005      3.812      0.000       0.009       0.028\n",
      "L1.e(dln_inv)    -0.2616      0.106     -2.478      0.013      -0.468      -0.055\n",
      "L1.e(dln_inc)     0.5085      0.630      0.807      0.420      -0.726       1.743\n",
      "L2.e(dln_inv)     0.0317      0.148      0.214      0.831      -0.259       0.322\n",
      "L2.e(dln_inc)     0.1905      0.476      0.400      0.689      -0.742       1.123\n",
      "                           Results for equation dln_inc                          \n",
      "=================================================================================\n",
      "                    coef    std err          z      P>|z|      [0.025      0.975]\n",
      "---------------------------------------------------------------------------------\n",
      "intercept         0.0207      0.002     13.018      0.000       0.018       0.024\n",
      "L1.e(dln_inv)     0.0476      0.042      1.141      0.254      -0.034       0.129\n",
      "L1.e(dln_inc)    -0.0709      0.140     -0.505      0.614      -0.346       0.204\n",
      "L2.e(dln_inv)     0.0181      0.042      0.427      0.670      -0.065       0.101\n",
      "L2.e(dln_inc)     0.1246      0.154      0.811      0.417      -0.177       0.426\n",
      "                             Error covariance matrix                              \n",
      "==================================================================================\n",
      "                     coef    std err          z      P>|z|      [0.025      0.975]\n",
      "----------------------------------------------------------------------------------\n",
      "sigma2.dln_inv     0.0020      0.000      7.359      0.000       0.001       0.003\n",
      "sigma2.dln_inc     0.0001   2.33e-05      5.829      0.000       9e-05       0.000\n",
      "==================================================================================\n",
      "\n",
      "Warnings:\n",
      "[1] Covariance matrix calculated using the outer product of gradients (complex-step).\n"
     ]
    }
   ],
   "source": [
    "mod = sm.tsa.VARMAX(endog[['dln_inv', 'dln_inc']], order=(0,2), error_cov_type='diagonal')\n",
    "res = mod.fit(maxiter=1000, disp=False)\n",
    "print(res.summary())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Caution: VARMA(p,q) specifications\n",
    "\n",
    "Although the model allows estimating VARMA(p,q) specifications, these models are not identified without additional restrictions on the representation matrices, which are not built-in. For this reason, it is recommended that the user proceed with error (and indeed a warning is issued when these models are specified). Nonetheless, they may in some circumstances provide useful information."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\ProgramData\\Anaconda3\\lib\\site-packages\\statsmodels\\tsa\\statespace\\varmax.py:159: EstimationWarning: Estimation of VARMA(p,q) models is not generically robust, due especially to identification issues.\n",
      "  EstimationWarning)\n",
      "C:\\ProgramData\\Anaconda3\\lib\\site-packages\\statsmodels\\tsa\\base\\tsa_model.py:165: ValueWarning: No frequency information was provided, so inferred frequency QS-OCT will be used.\n",
      "  % freq, ValueWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                             Statespace Model Results                             \n",
      "==================================================================================\n",
      "Dep. Variable:     ['dln_inv', 'dln_inc']   No. Observations:                   75\n",
      "Model:                         VARMA(1,1)   Log Likelihood                 354.284\n",
      "                              + intercept   AIC                           -682.567\n",
      "Date:                    Tue, 24 Dec 2019   BIC                           -652.440\n",
      "Time:                            15:05:56   HQIC                          -670.538\n",
      "Sample:                        04-01-1960                                         \n",
      "                             - 10-01-1978                                         \n",
      "Covariance Type:                      opg                                         \n",
      "===================================================================================\n",
      "Ljung-Box (Q):                68.75, 39.04   Jarque-Bera (JB):         10.79, 14.11\n",
      "Prob(Q):                        0.00, 0.51   Prob(JB):                   0.00, 0.00\n",
      "Heteroskedasticity (H):         0.43, 0.91   Skew:                      0.00, -0.46\n",
      "Prob(H) (two-sided):            0.04, 0.81   Kurtosis:                   4.86, 4.92\n",
      "                           Results for equation dln_inv                          \n",
      "=================================================================================\n",
      "                    coef    std err          z      P>|z|      [0.025      0.975]\n",
      "---------------------------------------------------------------------------------\n",
      "intercept         0.0110      0.068      0.162      0.872      -0.122       0.144\n",
      "L1.dln_inv       -0.0087      0.719     -0.012      0.990      -1.419       1.401\n",
      "L1.dln_inc        0.3609      2.844      0.127      0.899      -5.213       5.935\n",
      "L1.e(dln_inv)    -0.2508      0.731     -0.343      0.731      -1.683       1.182\n",
      "L1.e(dln_inc)     0.1263      3.086      0.041      0.967      -5.922       6.175\n",
      "                           Results for equation dln_inc                          \n",
      "=================================================================================\n",
      "                    coef    std err          z      P>|z|      [0.025      0.975]\n",
      "---------------------------------------------------------------------------------\n",
      "intercept         0.0165      0.029      0.580      0.562      -0.039       0.072\n",
      "L1.dln_inv       -0.0334      0.286     -0.117      0.907      -0.595       0.528\n",
      "L1.dln_inc        0.2318      1.160      0.200      0.842      -2.042       2.506\n",
      "L1.e(dln_inv)     0.0887      0.293      0.303      0.762      -0.485       0.663\n",
      "L1.e(dln_inc)    -0.2359      1.192     -0.198      0.843      -2.573       2.101\n",
      "                                  Error covariance matrix                                   \n",
      "============================================================================================\n",
      "                               coef    std err          z      P>|z|      [0.025      0.975]\n",
      "--------------------------------------------------------------------------------------------\n",
      "sqrt.var.dln_inv             0.0449      0.003     14.533      0.000       0.039       0.051\n",
      "sqrt.cov.dln_inv.dln_inc     0.0017      0.003      0.646      0.518      -0.003       0.007\n",
      "sqrt.var.dln_inc             0.0116      0.001     11.662      0.000       0.010       0.013\n",
      "============================================================================================\n",
      "\n",
      "Warnings:\n",
      "[1] Covariance matrix calculated using the outer product of gradients (complex-step).\n"
     ]
    }
   ],
   "source": [
    "mod = sm.tsa.VARMAX(endog[['dln_inv', 'dln_inc']], order=(1,1))\n",
    "res = mod.fit(maxiter=1000, disp=False)\n",
    "print(res.summary())"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
